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(V)  A one  dimensional  thermal  transport  code  has  been  employed  to  study  electron  energy  confinement  in 
tokamaks.  At  each  time  step  of  the  computations  each  point  in  the  profile  is  examined  to  determine  whether  it 
is  linearly  suble  to  either  the  dissipative  trapped  electron  mode  or  the  intenul  kink  and  tearing  modes.  If  the 
point  is  determined  to  be  unstable  an  anomalously  large  electron  thermal  conductivity  is  turned  on  in  the  vicinity 
of  that  point.  The  electron  temperature  profile  is  found  to  relax  to  a marginally  stable  state.  Our  code  predicts 
the  scaling  of  confinement  time  with  density  in  Alcator,  and  the  scaling  of  temperature  half-width  with  q(a)  in 
TFR. 
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A THEORY  OF  ELECTRON  ENERGY  CONFINEMENT  IN  TOKAMAKS 


I.  [ntroduction 

The  problem  of  electron  energy  flux  is  one  of  the  most  perplexing  problems  in  tokamak 
research.  The  most  recent  transport  codes,  have  focussed  on  the  six  regime  transport  model'. 
There,  depending  on  the  parameters,  the  plasma  is  subject  to  one  of  a number  of  instabilities. 
The  nonlinear  theory  of  each  instability  then  gives  rise  to  a transport  coefficient  in  the  ap- 
propriate regime.  Other  codes  have  focussed  on  the  so  called  pseudo  classical  model  for 
transport  ^ 

This  paper  examines  another  approach  to  the  problem  of  anomalous  transport  in 
tokamaks,  the  so  called  marginal  stability  approach  The  idea  is  that  as  soon  as  the  plasma 
becomes  unstable,  the  electron  transport  coefficient  jumps  abruptly  to  a large  value.  The 
resulting  enhanced  transport  forces  the  plasma  to  a stable  regime.  If  some  other  mechanism 
forces  the  plasma  towards  instability,  the  logical  steady  state  is  at  the  marginally  stable  point. 
The  basic  idea  behind  the  marginal  stability  approach  is  described  in  Ref.  8.  Recently,  we 
have  derived  relative  temperature  profiles  for  tokamaks^.  This  paper  is  an  extension  of  our 
earlier  work.  The  marginal  stability  approach  has  also  worked  well  in  describing  the  behavior 
of  resistive  shocks  the  behavior  of  plasma  in  electric  fields^’",  screw  pinches  and  also  the 
behavior  of  trapped  particles  in  the  magnetosphere'^.  Also,  there  is  apparantly  some  thermo- 
dynamic evidence  for  the  occurence  of  a marginally  stable  state  '^.  There  are  at  least  two  re- 
cent experimental  results,  the  qualitative  features  of  which  are  explained  very  naturally  be 
marginal  stability  theory.  The  first  is  the  result  of  heating  TFR  with  neutural  beams'^. 
Although  a large  amount  of  beam  power  was  deposited  into  the  electrons  the  electon  tempera- 
ture did  not  change  at  all  in  one  experiment,  and  changed  only  very  slightly  (less  than  10%)  in 
another.  This  result  is  certainly  consistent  with  the  idea  that  the  electron  temperature  profile 
is  resting  at  some  marginally  stable  value.  The  harder  it  is  pushed  into  the  unstable  regime, 
the  harder  the  plasma  fights  back.  The  result  is  that  the  temperature  profile  always  sits  near 
the  marginally  stable  profile  no  matter  how  hard  is  it  pushed. 

The  second  is  the  results  of  heat  pulse  propagation  studies  in  ORMAK'*.  A pulse  of  en- 
ergy is  dumped  from  the  center  of  the  plasma  rja  = 0 (a  — limiter  radius)  to  about  rja  = .25. 
The  propagation  of  the  heat  pulse  to  rja  = .75  is  studied  by  X-rays.  The  remarkable  fact  is 
that  the  time  for  pulse  propagation  outward  is  considerably  shorter  than  the  energy  contain- 
ment time.  This  is  so  even  though  the  relative  temperature  perturbation  is  quite  small.  This 
too  is  consistent  with  the  plasma  temperature  remaining  at  a marginally  stable  profile;  any  in- 
crease in  temperature  drives  the  plasma  to  an  unstable  state,  and  is  therefore  rapidly  expelled. 

To  examine  whether  the  concept  of  marginal  stability  is  viable  for  a tokamak  plasma,  we 
have  put  together  a small  transport  code  which  solves  for  electron  and  ion  temperatures.  A 
large  electron  thermal  conductivity  is  turned  on  whenever  the  plasma  becomes  linearly  un- 
stable. Our  basic  assumption  is  that  in  a real  Tokamak  the  anaomalous  electron  thermal  con- 
ductivity that  would  result  if  the  plasma  were  forced  into  an  unstable  state  is  large  enough  to 
maintain  the  temperature  near  a linearly  marginaly  stable  profile.  Thus,  as  far  as  determining 
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the  profile  is  concerned  all  that  is  necessary  to  know  is  the  marginal  stability  condition  (the 
conditions  under  which  the  transport  should  become  large).  The  magnitude  of  the  transport  in 
an  unstable  plasma  affects  only  the  rate  at  which  the  plasma  is  driven  to  the  marginally  stable 
state. 


Our  approach  is  similar  to  a quasilinear  theory  of  transport  with  the  additional  assumption 
that  the  characteristic  time  over  which  the  instability  can  develoo  and  field  amplitudes  vary  is 
short  compared  to  the  time  scale  of  the  gross  energy  confinement  in  the  plasma.  Thus  if  the 
field  amplitudes  are  quasisteady  indicating  saturation  of  the  instability  the  linear  growth  rate 
must  be  nearly  zero.  This  code  however  does  not  model  as  many  physical  effects  as  does 
others'.  Specifically  we  do  not  solve  for  the  density  profile,  the  impurity  transport  or  the  neu- 
tural  transport.  Also,  we  assume  that  the  current  is  always  its  steady  state  value 
U — <r  so  that  the  code  does  not  solve  for  skin  effect  either.  Our  goal  is  only  to  see 

whether  the  marginal  stability  approach  is  a viable  concept.  If  it  is,  a more  complete  code  can 
be  extended  or  modified  to  incorporate  it. 

In  order  to  utilize  this  approach,  it  is  necessary  to  carefully  examine  the  stable  and  un- 
stable regimes  of  parameter  space.  We  find  that  the  instability  does  not  depend  only  on  col- 
lisionality.  The  density  gradient,  temperature  gradient  and  shear  all  play  important  roles.  One 
particularly  important  effect  is  that  the  temperature  gradient  can  be  strongly  stabilizing  at  low 
enough  temperatures,  so  that  the  outer  region  of  the  plasma  is  generally  classical.  At  higher 
temperatures  only  the  shear  is  stabilizing,  and  the  electron  temperature  picks  out  the  profile  at 
which  shear  stabilization  just  balances  the  temperature  and  density  gradient  induced  growth 
rates. 


The  results  of  our  code  are  sufficiently  encouraging  that  we  conclude  that  the  marginal 
stability  approach  is  viable.  Its  theoretical  foundation  is  very  simple  and  understandable;  it 
does  not  rely  at  all  on  complex,  speculative  nonlinear  theories  of  many  different  instabilities. 
The  only  uncertainties  involved  are  those  inherent  in  linear  theory. 

As  far  as  confinement  times  are  concerned,  our  predictions  are  generally  too  high  by  a 
factor  of  from  two  to  five.  However  temperature  profiles  and  scaling  laws  show  much  better 
agreement  with  experiment.  Two  scaling  laws  which  our  code  predicts  are  the  increase  of 
confinement  time  with  density  in  Alcator,  and  the  fact  that  the  temperature  half  width  and  the 
radius  of  the  9 — 1 surface  depend  only  on  a single  parameter,  as  measured  in  TFR. 
Furthermore,  our  calculations  of  fluctuation  levels  also  agree  resonably  with  those  measured  in 
ATC. 


II.  The  Code 
A — Basic  Fluid  Code 


In  order  to  test  whether  marginal  stability  is  a reasonable  approach  to  the  problem  of 
anomalous  electron  energy  transport  in  tokamaks,  we  have  constructed  a simple  computer 
code  The  equations  which  are  solved  numerically  are 
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where  T",.,,,  are  the  electron  (ion  temperature).  7 is  the  current  density,  £ is  the  toroidal  elec- 
tric field,  assumed  to  be  independent  of  radius,  ir  is  the  electrical  conductivity  given  by 
Rutherford ', 
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is  the  electron  ion  momentum  exchange  collision  time 
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where  is  the  effective  charge  of  the  plasma.  Also,  and  are  the  neoclassical  expres- 
sion for  electron  and  ion  thermal  conductivity.  Here  was  adopt  the  expression  given  by  Ruth- 
erford, Duchs*’^, 
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Above,  (p,„  ) is  the  electron  (ion)  larmor  radius  in  the  poloidal  field,  t,  is  the  ion-ion  colli- 
sion time 


Ti  — (3/w/^^  ) (4Vff  w'*ln  .\,  ) '.  (5) 

The  quantities  i',  and  are  the  ratio  of  effective  ion  or  electron  collision  frequency  (r« ) 
to  the  ion  or  electron  bounce  frequency 

«'/  “ ^ef^  ^TJm)  (6a) 

v'  - Z„  ^1/2  (7-y^)  1/2^  (^J,) 

and  all  other  notation  is  standard.  Notice  that  the  fluid  equations  solved  have  no  magnetic 
diffusion,  but  rather  assume  the  steady  state  value  for  current.  Thus,  the  code  does  not 
correctly  predict  the  approach  to  equilibrium.  However  it  should  predict  the  correct  steady 
state.  Finally,  is  the  anomalous  electron  thermal  conductivity.  Consistent  with  our  as- 
sumption of  marginal  stability,  we  take  for  a function  which  rises  very  suddenly  at  the 
onset  of  instability. 


There  are  two  instabilities  we  consider  here.  The  first  are  internal  kink  and  tearing 
modes  which  occur  at  the  center  of  the  plasma.  These  are  observed  to  have  an  w — « — 1 
structure  and  are  detected  by  saw-tooth  oscillations  on  the  soft  X-ray  signals^®’^'.  To  simulate 
the  effect  of  these  modes  we  assume  jumps  to  a large  value  whenever  the  safety  factor 
q - rB/RB„  is  less  than  unity.  This  is  roughly  the  same  approach  taken  in  the  Princeton  nu- 
merical simulations'. 


The  second  instability,  which  we  concentrate  on  here  is  the  dissipative  trapped  electron 
mode^^.  This  will  be  discussed  more  fully  in  a later  subsection.  For  our  purpose  now,  we  can 
regard  the  plasma  as  unstable  to  this  instability  if  the  ratio  of  density  gradient  scale  length  to 
shear  scale  length,  H,  is  less  than  some  critical  value  (4^^.  The  form  for  the  anomalous  thermal 
conduction  which  we  choose  is 
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1 + <7 


1 + («/«,,)  ■" 


where  D/j  is  the  Bohm  diffusion  coefllcient 

and  n,  /;>>!.  We  have  called  the  functions  in  the  parentheses  thermostat  functions  be- 
cause they  jump  sharply  from  zero  to  one  as  q and  go  from  slightly  greater  than  one  to 

slightly  less  than  one. 

Many  of  the  quantities  in  Eqs.  (1)  — (8)  depend  on  density.  Therefore  to  run  the  code, 
it  is  necessary  to  specify  a density  profile,  in  a later  subsection,  we  discuss  how  the  code  could 
be  extended  to  solve  for  the  density  profiles  self  consistently.  For  all  of  the  work  presented 
here  however,  density  profiles  are  specified  arbitrarially. 

The  remaining  things  to  specify  are  boundary  conditions  for  temperature  and  current. 
We  choose  an  electron  and  ion  temperature  of  10  ev,  at  the  limiter  and  consider  the  total 
current  to  be  specified.  This  is  equivalent  to  choosing  a particular  value  of  q at  the  limiter. 

B — The  Role  of  Nonlinear  Theory 

It  is  worth  reiterating  here  that  there  is  a very  fundamental  difference  in  the  approach 
taken  here  and  the  approach  taken  by  others.  Generally  plasma  is  assumed  to  exist  in  a linear- 
ly unstable  state.  The  transport  coefficient  is  determined  by  the  fluctuation  level  which  is  lim- 
ited by  some  local  nonlinear  effect,  for  instance  resonance  broadening  mode  coupling 
etc.  Therefore,  to  get  transport  coefficients,  it  is  essential  to  have  a nonlinear  theory  of  all 
relevant  instabilities.  In  our  approach,  we  assume  the  most  effective  stabilization  mechanism  is 
the  relaxation  of  the  plasma  to  a linearly  stable  state.  If  some  mechanism  forces  the  plasma  to- 
wards instability,  it  is  reasonable  to  assume  that  the  plasma  sits  at  marginal  stability.  For  an 
Ohmically  heated  tokamak,  this  mechanism  is  the  channeling  of  the  current  into  the  hotter  re- 
gions in  the  center.  This  then  increases  the  electron  temperature  gradient  (since  ir  — a)  is 
held  fixed  at  lOev)).  However  a steep  electron  temperature  gradient  drives  trapped  particle  in- 
stabilities which  then  causes  anomalous  thermal  conduction.  Thus  a dynamic  balance  should 
be  possible,  with  the  tokamak  temperature  profile  sitting  at  marginal  stability. 

Needless  to  say,  nonlinear  theory  of  the  relevant  instabilities,  plays  a much  less  important 
role  in  this  type  of  plasma,  than  it  does  if  the  plasma  is  assumed  to  be  everywhere  linearly  un- 
stable. This  is  not  to  say  that  nonlinear  theory  plays  no  role  at  all.  Indeed,  the  proper  evalua- 
tion of  the  thermostat  functions,  defined  in  the  last  section,  can  only  be  made  through  non- 
linear theory.  It  is  interesting  to  note  however  that  the  most  recent  studies  of  stabilization  of 
drift  waves  and  trapped  particle  instabilities  by  resonance  broadening  show  no  nonlinear  sa- 
turation unless  y/<u  < (k||V,/o>)^.  That  is,  no  saturation  unless  the  mode  is  only  slightly 
above  marginal  stability. 


C — The  Instability 

It  remains  to  determine  H This  can  be  done  by  examining  the  dispersion  relation  for 
trapped  electron  instability  The  theory  of  this  instability  is  continually  being  refined, 

and  it  would  be  difficult  to  write  an  expression  for  growth  rate  which  would  be  universally  ac- 
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cepled  at  this  lime.  The  disperion  relation  we  use  neglects  the  coupling  between  poloidal  and 
radial  structure  so  it  does  not  properly  treat  the  two  dimensional  mode  structure^'.  However  it 
probably  does  reasonably  approximate  all  of  the  significant  physical  effects. 


The  dispersion  relation  we  use  is^^ 

* - .S(b)lr  II  -7,,*^  G(b)ll  -1  -r 
to 
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In  Eq.  (9)  above  t ~ TJT  , £ “hi./hi,  -v^/(v^t^€)  . toj  -hi.«„vVv^, 
« — HR,  Sib)  - !„  ib)e  “*  b — k^ft} . v,  is  the  density  gradient  scale  length  divided  by  the 
ion  temperature  gradient  scale  length,  Gib)  — 1 — /,  ib^)//g  ib^),  UT  denotes  an  integral 
over  untrapped  regions  of  electron  velocity  space  (v||/v|  > V2«  and  also  ni^fv)  > 
ir^/f  ) (Vy/v)^  where  toi,iv)  •«  s/2e  v/Rq),  and  T" denotes  the  integral  over  the  trapped  region 
of  electron  velocity  space.  The  quantity  hi.^  — to*  (1  + Tig  (1/2  - 3/2))  where  HI.  is 

the  electron  diamagnetic  drift  velocity  and  is  the  ion  collisional  damping  rate^^ 
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To  derive  Eq.  (9),  we  have  split  the  electron  velocity  space  up  into  trapped  and  untrapped  re- 
gions. The  to.tg  (v/v^)^  term  in  the  denominator  of  the  trapped  term  is  the  effect  of  the  drift 
resonances.  Here  untrapped  term,  we  have  assumed  A||  — \/Rq.  To  cal- 

culate the  effect  of  shear  induced  damping,  we  have  neglected  the  poloidal  mode  structure  and 
have  performed  a Pearlstein-Berk  calculation  of  the  shear  damping  assuming  A,.p,  < < 1 
and  also  assuming  that  a mode  unstable  on  one  mode  rational  surface  does  not  couple  to  waves 
on  other  rational  surfaces.  Unfortunately,  this  approximation  is  not  valid  for  A„p,  larger  than 
about  0.2.  For  larger  Ap,,  the  calculation  of  shear  induced  damping  is  extremely  complicated, 
involving  the  full  two  dimensional  nature  of  the  mode  structure^'.  Preliminary  indications  are 
that  shear  is  somewhat  less  effective  a stabilization  mechanism  than  is  described  in  Eq. 

(9).  We  have  arbitrarily  reduced  the  effect  of  shear  stabilization  in  our  code  in  certain  cir- 
cumstances by  multiplying  the  last  term  in  Eq.  (9)  by  a factor  less  than  unity  in  order  to  assess 
the  consequences  of  this  possibility. 


We  will  now  discuss  qualitatively  the  effects  of  various  terms  in  the  dispersion  relation. 

1)  In  both  the  T and  UT  terms,  the  effect  of  a density  gradient  is  destabilizing. 

2)  Shear  is  stabilizing  for  all  A and  all  temperature. 

3)  An  electron  temperature  gradient  is  stabilizing  in  the  UT  term. 

4)  An  electron  temperature  gradient  is  destabilizing  in  the  T term. 
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S)  Ion-ion  collisions  are  always  stabilizing  and  can  be  an  important  stabilizing  effect  at 
low  enough  temperature. 

From  these  facts,  some  rough  qualitative  features  of  the  instability  can  be  discerned.  As 
the  plasma  temperature  is  lowered,  there  are  fewer  and  fewer  trappped  particles  to  drive  insta- 
bility. Therefore  as  the  temperature  is  lowered,  the  electron  temperature  gradient  will  at  some 
point  change  from  a destabilizing  to  a stabilizing  effect.  Since  temperature  profiles  in  tokamaks 
are  generally  steeper  than  density  profiles,  we  expect  the  outer,  cooler  regions  of  the  tokamak 
to  be  stable.  As  a very  rough  rule  of  thumb,  we  find  that  the  electron  temperature  gradient 
changes  from  a destabilizing  to  a stabilizing  effect  at  approximately  the  transition  from  the  pla- 
teau regime  to  the  Pfirsch-Schliter  regime. 

The  next  problem  is  to  analyize  the  dispersion  relation  to  find  the  critical  shear  for  stabil- 
ity. To  make  any  progress,  it  is  necessary  to  make  a few  approximations,  principally  to  shrink 
down  the  dimensionality  of  the  parameter  space.  We  wish  to  determine  from  Eq.  (9)  the  value 
of  W,  the  shear,  necessary  to  stabilize  all  modes,  that  is  we  require  /m(a)(/cp,  )l  < 0 for  all  fcp,. 
To  accomplish  this  we  calculate  the  value  of  B necessary  to  stabilize  the  mode  for  each  partic- 
ular itp,,  call  it  » (kpj),  and  select  B^^  to  be  the  largest  of  these  values.  Thus,  we  view  the  real 
and  imaginary  parts  of  Eq^^  (9)  as  two  equations  determing  the  two  unknowns 
to^  (oi  “ tu,  -F  iy,  y “ 0)  and  # (ftp,). 

The  real  part  of  Eq.  (9)  does  not  involve  B and  hence  may  be  solved  to  determine  <«>,. 
Once  £u,  is  known  for  each  kp,  it  is  a simple  matter  to  substitute  tu,  (Arp,)  in  the  imaginary  part 
of  Eq.  (9)  to  determine  B(ilcp,).  Here  we  point  out  that  <i>,.  is  the  real  frequency  for  a margi- 
nally stable  mode  and  not  necessarily  the  real  part  of  the  frequency  for  a damped  mode.  In 
our  computations  we  assumed  that  the  nonadiabatic  electron  response  was  negligible  insofar  as 
determining  to,.  In  this  case  the  equation  for  to,  is  given  by, 

S(6)(T-Fni  -yi.b^GW]]  - 1 - r - 0,  (11) 

where  r — and  ? — to  Vto.  Equation  (11)  was  further  simplied  by  assuming  t - 1.5, 
and  Tj,  — 1.0.  This  allowed  us  to  calculate  to/to*  as  a function  of  Arp,  from  the  outset  and, 
hence,  it  was  not  necessary  to  solve  Eq.  (11)  at  each  point  in  the  profile  or  each  time  step. 

The  imaginary  part  of  Eq.  (9)  can  be  written  as  follows, 

B I (b)  -/,  (o/,  6„.  e,  4,  ?)  + t)^/2  (o^*,  e„,  6.  4,  (?)  - 8„,  (12) 

where 

Z “ 4 rSib)  1 It  4 (1  — tj;  ) 1 G(.b) 

-h  17/4  (1  ~ (1  — 28^)  C(8)  ] '^^1  It  -F  4 tl  -Ft?,  — i7,7>^  C7(6)  ] j (13) 

and  /|  and  /2  are  functions  involving  integrals  over  velocity  space  that  result  from  the  nonadi- 
abatic electron  response. 

The  first  term  in  Eq.  (12)  is  the  density  gradient  driving  term.  The  function  /,  is  always 
greater  than  zero,  indicating  that  the  density  gradient  is  always  destabilizing.  The  second  term 
is  the  temperature  gradient  driving  term.  The  function  can  be  either  negative  or  positive 
depending  on  whether  the  temperature  gradient  is  stabilizing  (Pfirsch  Schliiter  regime  accord- 
ing to  our  rule  of  thumb),  or  destabilizing  (plateau  and  banana  regime).  The  third  term  is  the 
ion-ion  collisional  damping.  We  use  polynomial  approximations  for  J\ , and  /2,  which  are  dis- 
cussed fully  in  Appendix  A.  Then,  from  Eq.  (13),  we  can  solve  for  the  critical  shear,  B„,  for 
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each  kit,.  Once  Wj,  has  been  determined  we  have  everything  needed  to  advance  the  fluid 
quantites  a time  step. 


In  our  computation  Eqs.  (11)  and  (12)  were  solved  for  five  values  of  kf),.  These  values 
» 0.25,  0.5,  0.75,  1.0,  1.25.  In  order  to  determine  the  fluctuation  level, 
eth/T,.  - f>n/ng,  we  simply  see  what  fluctuation  level  is  required  to  produce  the  necessary  tran- 
sport. The  quasi  linear  expression  for  the  hermal  conductivity  is 


T — 

r 3^ 


20  .2  l<f>(A)l^c^  y/2^  , 

m ^ • n. 


Bi  T, 


j —3) 

tit  + iVfj  - utj 


(14) 


where  the  triangular  brackets  mean  an  average  over  the  magnitude  of  the  particle  velocity. 
Thus  from  the  known  value  of  Kg„  and  Ap,,  a fluctuation  level  may  be  calculated. 


The  Density  Response 

Although  we  only  solve  for  the  temperature  profile  with  an  assumed  electron  density 
profile,  it  is  possible  to  solve  for  the  density  profile  also.  The  resulting  calculations  however 
would  be  much  more  involved  for  at  least  two'  reasons.  First  of  all,  using  a specified  density 
profile  for  all  time  greatly  simplified  the  numerics  involved  in  the  solution  of  the  dispersion  re- 
lation. Secondly,  to  solve  for  the  density,  one  must  include  all  sources  of  electrons.  This 
means  solving  the  coupled  system  of  ion,  neutral,  and  impurity  transport,  and  also  utilizing  re- 
liable models  for  ion,  neutral,  and  impurity  recyclying  at  the  walls. 


However,  in  principle  at  least,  it  is  not  difficult  to  solve  for  the  density  the  same  way  as 
we  solve  for  the  temperature.  The  density  equation  is 


8w 

dt 


(15) 


where  is  the  neoclassical  diffusion  coefficient,  5 is  a term  for  all  sources  and  sinks  of  elec- 
tron density  from  Ware  Pinch,  ionization  and  recombination  of  impurities,  and  ionization  and 
recombination  of  neutrals.  Approximate  expressions  for  S can  be  found  in  the  Princeton  work. 
Finally  is  the  anomalous  Diffusion  coefficient 


1 + <l 


■r-  + 


(0/e^,) 

1 -t-  (e/e„)  ■ 


(16) 


In  Eq.  (16)  above  the  factors  a and  a'  result  from  the  fact  that  the  thermal  conduction  and 
diffusion  coefficients  are  not  the  same.  The  quasi-linear  theory  of  the  trapped  electron 
instability^^  gives  a value  for  a.  This  quantity  usually  has  a value  less  than,  but  of  order  unity, 
typically  a is  about  one  half  or  one  third.  (Actually,  for  the  trapped  particle  instability  the 
temperature  gradients  can  drive  particle  fluxes  and  the  density  gradients  can  drive  energy 
fluxes.  These  effects  could  also  be  included  in  Eqs.  (1)  and  (15).)  An  expression  for  n' would 
have  to  come  from  the  theory  of  the  MHD  instabilities  occurring  near  the  9 “1  surface  in 
the  center  of  the  tokamak. 


To  summarize,  it  is  possible  to  extend  the  work  presented  here  to  include  the  density 
response  as  well  as  the  temperature  response.  Whenever  the  plasma  is  unstable,  anomalous 
diffusion  is  turned  on  and  the  density  profile  will  quickly  relax.  The  anomalous  diffusion 
coefficient  is  then  related  to  the  anomalous  thermal  conduction  coefficient  by  the  quasi-linear 
theory  of  the  relevant  instability. 
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The  Dependence  of  the  Temperature  Profile  on 

the  Density  Profile  and  Exponent  of  the  Thermostat  Function 

In  this  subsection  and  the  next  few  subsections  we  describe  a series  of  tests  of  the  code. 
For  all  tests,  we  have  used  PLT  parameters  in  experimental  runs  recently  reported These 
are  R -135  cm,  a - 45  cm,  qia)  -4.6,  - 35  KG,  - 4 and  n{r  -0)  -5 

X 10  cm  The  first  run  used  exponents  n - 8,  //'  — 20  and  a density  profile  n ^ — 

0. 85  (r/a)^).  The  steady  state  electron  temperature  is  shown  .n  the  dashed  curve  of  Fig.  1. 
Clearly,  there  is  a very  large  unstable  region  of  the  plasma,  extending  nearly  to  the  edge.  The 
density  profile  we  have  chosen  is  so  unstable  that  no  amount  of  shear  consistent  with  q(0)  — 

1,  <((a)  — 4.6  can  stabilize  it.  Thus  the  temperature  gradient  flattens  out  over  a very  wide  re- 

gion and  the  heat  flux  is  characterized  by  Bohm  diflusion.  Finally,  near  enough  to  the  edge, 
the  plasma  is  linearly  stable  and  the  diflusion  is  classicial,  but  with  a very  large  temperature 
gradient.  This  result  and  the  physical  model  used  to  derive  it  are  obviously  unphysical  for  a 
number  of  reasons.  First  of  all,  an  electron  temperature  of  2 KEV  so  near  the  limiter  is  almost 
certainly  unphysical.  Second,  it  is  not  clear  how  accurate  the  numerics  are  for  such  a profile; 
there  are  only  two  grid  cells  in  the  steep  temperature  drop  region  at  the  edge  (although  the 
final  steady  state  solution  does  remain  constant  in  time).  Third,  and  most  important  from  our 
point  of  view,  the  density  profile  is  assumed  frozen  into  a very  unstable  shape.  The  instability 
will  in  fact  relax  this  density  profile,  as  discussed  in  the  previous  subsection.  Since  our  code 
does  not  allow  the  density  profile  to  relax,  we  have  tried  another  run  with  the  same  density 
profile,  but  where  the  density  profile  cannot  itself  drive  instability.  That  is  we  have  set  /j  — 0, 
but  with  all  density  gradient  effects  (for  instance  drift  resonances)  included  in  In  a very 
crude  way,  this  mocks  up  the  fact  that  the  density  gradient  has  relaxed  and  can  no  longer  drive 
the  instability.  The  resulting  electron  temperature  profile  is  shown  as  the  dash  dot  curve  m 
Fig.  (1).  Notice  that  the  temperature  profile  now  looks  much  more  resonable  and  physicv* 

The  fact  that  the  temperature  peak  is  somewhat  higher  is  also  resonable  since  the  profile  is 

now  more  stable. 

As  a final  example  of  the  effect  of  density  profile  on  temperature,  we  ran  the  code  with  a 

more  stable  density  profile  n — ng(l  —0.5  (r/a)^).  The  result  is  shown  in  Fig.  1 as  the  solid 

curve.  Since  this  seems  to  be  the  most  reasonable  temperature  profile,  in  all  our  subsequent 
runs,  we  used  this  density  profile  and  retained  /2  in  the  dispersion  relation. 

The  next  thing  to  examine  is  the  dependence  of  the  result  on  the  exponent  of  the  ther- 
mostat function.  Clearly,  as  /i  — <»  the  temperature  profile  will  be  independent  of  n.  Actually 
the  shape  of  the  thermostat  function  depends  upon  the  nonlinear  theory  of  the  instability.  It 
may  be  that  the  effect  of  nonlinear  theory  (as  manifest  in  the  shape  of  the  thermostat  func- 
tion) plays  an  important  role.  For  our  purposes  , we  neglected  this  aspect  of  the  anomalous 
transport  and  chose  an  n large  enough  that  the  temperature  profile  is  independent  of  n.  The 
steady  state  temperature  profiles  for  3 different  values  of  n,  n —4,  16  and  32  are  shown  in  Fig. 

2.  (The  temperature  profile  for  w — 8 is  shown  in  Fig.  1).  Notice  that  the  result  for  « — 4 
gives  a somewhat  lower  temperature  than  the  results  for  n — 16  and  32.  The  reason  is  that  in 
both  cases,  the  thermal  conduction  is  a very  small  fraction  of  the  Bohm  value.  However  for  (-) 
< <■>...  the  thermostat  function  for  n — 4 is  larger  than  that  for  n —32,  so  the  thermal  conduc- 
tion is  greater  and  the  temperature  should  be  smaller.  Henceforth,  we  retained  the  choice  n — 
32  in  all  runs. 

As  a final  check  on  the  results  of  the  code,  we  have  taken  the  solid  temperature  profile  of 
Fig.  1 and  have  numerically  solved  the  full  linear  dispersion  relation,  Eq.  (9),  at  4 points  on 
this  profile.  The  numerical  solution  of  the  dispersion  relation  confirmed  that  the  profile  is  mar- 
ginally stable.  At  the  points  in  the  profile  where  the  M thermostat  function  is  zero,  the  numer- 
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ical  soululion  of  Eq.  (9)  showed  no  instability.  At  the  points  where  the  <-)  thermostat  function 
is  greater  than  zero,  the  numerical  solution  of  Eq.  (9)  showed  a very  slighly  growing  or  damped 
roots.  The  plasma  therefore  is  very  close  to  marginal  stability. 

Classical  Versus  Neo-Classical  Resistivity 

One  possible  question  is  whether  to  use  the  classical  expression  for  conductivity  or  the 
neo-classical  expression 

<r,vc  -"-(I  - (1.9(r/f?)'''2  -rlR)l(\  -E  *<;)).  (17) 

In  Figs.  3a  and  b are  shown  plots  of  electron  temperature  profile  and  q(r)  profile  for  the  use  of 
classical  and  nonclassical  resisitivty.  Notice  that  the  temperature  profiles  for  the  two  are  quite 
similar.  However  the  q profiles  are  strikingly  different.  For  the  case  of  classical  resisitivity  q is 
everywhere  above  one.  For  the  case  of  non-classicial  resisitivty,  q sinks  to  values  considerably 
smaller  than  one. 

The  reason  is  that  for  the  case  of  classicial  resisitivity,  the  resistivity  depends  only  on 
temperature,  so  that  flattening  the  temperature  profile  flattens  the  current  profile.  For  neo- 
classicial  resisitivity  however  the  resisitivity  depends  on  both  temperature  and  radius.  Thus 
flattening  the  electron  temperature  profile  (which  is  all  that  anomalous  thermal  conduction  can 
do)  does  not  flatten  the  current  profile.  Thus  the  instability  is  not  reacting  back  directly  on  the 
driving  mechanism.  (We  are  not  sure  how  to  interpret  this  result.  If  there  is  a strongly  exict- 
ed  internal  kink  or  tearing  mode  in  the  center,  the  magnetic  field  lines  will  almqst  certainly  not 
stay  on  smooth  circular  magnetic  surfaces,  but  rather  will  more  likely  wander  ergodigally 
throughout  the  region  q < 1.  It  is  not  at  all  clear  to  us  whether  neo-classical  expressions  for 
resisitivity  are  valid  in  this  case.  Possibly,  the  neoclasical  expression  for  resisitivity  should  be 
averaged  in  some  way  over  the  q < 1 region.  We  have  decided  to  avoid  this  problem  by  using 
the  ordinary  classical  expression  for  conductivity  (Eq.  (2))  so  that  anomalous  thermal  conduc- 
tion alone  flattens  out  both  the  temperature  and  current  profile.  However  this  question  may  be 
worth  additional  study.) 

Dependence  on  Shear  Damping 

Recently  there  has  been  some  speculation^'"’^"’*  to  the  effect  that  because  of  coupling  of 
different  mode  rational  surfaces,  shear  is  not  as  effective  a damping  mechanism  as  indicated  in 
the  linear  dispersion  relation  Eq.  (9).  We  have  run  the  code  with  the  shear  induced  damping 
reduced  by  a factor  of  two  (W  — »/2  in  Eq.  (9)).  The  plot  of  electron  temperature  is  shown  in 
Fig.  4.  Notice  that  the  central  temperature  drops  from  about  3.2  KEV  to  about  2.4  KEV.  In 
all  our  subsequent  runs,  we  retain  this  choice,  F)  — 0/2,  in  Eq.  (9). 

The  Energy  Confinement  Time 

As  will  be  shown  in  the  next  section,  jur  code  generally  gives  good  agreement  with  ex- 
periment on  central  temperature  and  temperature  profile,  but  not  nearly  as  good  agreement  on 
the  energy  confinement  time  The  reason  is  the  extreme  sensitivity  of  on  relevant 
parameters  for  ohmically  heated  tokamaks.  The  energy  confinement  time  is  defined  as  the  en- 
ergy content  divided  by  the  power  into  the  tokamak  IV  (current  times  voltage), 

f cl\ll2  n(T,  -I-  T,) 

r,  ^ (18) 

Let  us  call  the  central  temperature  T„  and  the  temper, ’nre  half  width  ^r(T).  Then  the  energy 
content  scales  as  Tg  Ar^fT"). 
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In  our  code,  only  the  current  is  specified,  the  voltage  must  be  calculated  from  the  resisi- 
tivity  profile.  The  conductivity  scales  as  and  the  total  resistance  scales  as  \r  (T).  Thus 
the  voltage  scales  as  lr~^(T).  Hence  the  theoretical  energy  confinement  time  scales 

as 

T^(Th)-  IrUT).  (19) 

Thus  a fifteen  percent  error  in  Tg  and  \r(T)  implies  an  error  of  more  than  a factor  of  two  in 
T p.  This  sensitivity  of  Tp  is  characteristic  not  only  of  our  code,  but  is  also  characteristic  of  any 
c^e  which  attempts  to  calculate  temperature  profiles  of  ohmically  heated  tokamaks  with 
specified  current. 

It  might  be  argued  that  such  a sensitivity  is  also  true  of  experimental  determinations  of 
Tp.  While  there  are  experimental  uncertainties  in  Tg  and  lr{T),  an  experimentalist  at  least 
has  an  accurate  measure  of  the  voltage.  Thus  experimentally,  the  confinement  time  scales  as 

Tp(EX)  ' Tg  Ir^iT).  (20) 

which  exhibits  much  less  sensitivity  to  the  temperature  and  temperature  profile  than  the  scal- 
ing shown  in  Eq.  (19).  Thus  we  conclude  that  relatively  inaccurate  predictions  of  for  ohmi- 
cally heated  tokamaks  is  something  transport  codes  must  live  with. 

Summary 

To  summarize,  we  assume  in  all  subsequent  runs  the  density  profile  n — «„(!  — 0.5 
(r/a)^),  classical  rather  than  neo-classical  resisitivty,  //'  — 20,  n — 32,  and  that  the  shear  damp- 
ing is  reduced  by  a factor  of  two.  For  PLT  parameters  /t  — 135  cm,  <7—45  cm,  ^(a)  —4.6,  it 
— 35  KG,  — 4,  Hg  — 4.2  x 10  and  a hydrogen  plasma  is  assumed,  the  prediction  of  our 
code  is  shown  in  Fig.  4.  Plotted  are  electron  and  ion  temperature,  q(r)  and  the  radial  depen- 
dence of  the  (/and  H thermostat  functions.  The  confinement  time  (electron  plus  ion)  is  200 
miliseconds.  Clearly  the  inner  region  is  marginally  stable  to  the  internal  kink  mode,  the  region 
between  Ha  —0.25  and  r —0.9  is  marginally  stable  to  the  trapped  electron  instability,  and  the 
outer  most  region  is  classical. 

III.  Results  for  TFR 

One  very  interesting  experimental  result  from  TFR^'  is  that  the  temperature  half  width 
Xr{T)  depends  on  only  a single  parameter,  q{a),  even  though  two  parameters,  fiand  / are  in- 
dependently varied.  This  result  comes  out  of  our  original  simplified  marginal  stability  theory^. 
There  the  driving  term  for  instability  depended  only  on  t),  and  r,  while  the  shear  stabilization 
term  depended  only  on  W and  r.  Only  the  relative  temperature  profile  was  solved.  However,  in 
the  simplified  theory,  the  only  parameter  characterzing  the  relative  temperature  profile  was 
q(a).  Thus  IriT)  depended  only  on  qia)',  the  larger  q{a),  the  more  peaked  the  temperature 
profile.  Although  we  now  use  a much  more  complicated  dispersion  relation,  the  same  basic 
result  still  holds  true.  In  Fig.  5 is  shown  (solid  line)  the  predictions  by  our  code  for  SriT)  as  a 
function  of  q (a). 

Here  we  have  assumed  Z^/  — 3 and  “ 6 x 10 and  a helium  plasma  was  assumed 
since  most  of  the  experiments  in  Ref.  21  used  a helium  plasma.  Various  dots  show  the  calcu- 
lated points.  Although  there  is  some  scatter  to  the  points,  ArlD  basically  does  depend  on  the 
single  parameter  q(a).  The  dashed  line  shows  the  experimental  result  from  TFR.  Note  that 
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our  code  predicts  somewhat  wider  temperature  half  widths  than  are  measured.  A possible  rea- 
son for  this  discrepancy  is  the  effect  of  radiation.  This  will  be  discussed  qualtatively  shortly. 

Another  interesting  result  from  the  TFR  experiment  is  that  the  position  of  the  </  “ 1 sur- 
face, as  defined  by  the  radial  position  of  the  mode  in  the  saw  tooth  oscillations  in  the  soft 
X/ray  signal,  also  depends  only  on  the  single  parameter  qia).  In  Fig.  S,  the  dash  dot  line 
shows  that  the  predicted  value  for  the  q “ 1 surface,  as  defined  by  the  radial  position  where 
the  q thermostat  function  drops  to  half  its  maximum  value,  also  depends  principally  on  this 
single  parameter. 

In  Fig.  6a  and  b are  predictions  of  the  code  for  electron  temperature  compared  to 
experiment*’,  for  three  discharges.  The  predicted  energy  confinement  times  (electron  plus 
ion)  are  7 p — 65,  96,  108  ms.  for  the  3 cases  B — 25,  40  and  50  Kg.  These  are  larger  than  the 
experimental  values  by  about  a factor  of  5.  Notice  however  that  the  code  does  give  good 
agreement  with  experiment  as  regards  central  temperature.  In  Fig.  7 are  plotted  the  radial 
dependence  of  T^„  Ti,  q,  for  the  run  with  B — 50  Kg,  / “ 140  KA. 

Other  quantities  of  interest  are  the  of  the  fluctuation  and  predicted  values  of  ed>/r^  as 
a function  of  radius.  These  are  shown  in  Fig.  8 for  a discharge  with  60  KG,  300  KA.  The 
fluctuation  level  seems  to  be  comparable  to,  but  somewhat  larger  than  recently  reported^® 
Actually  the  values  for  eiblT  given  in  Ref.  38  should  be  doubled'*®.  Notice  that  the  discon- 
tinuities in  etb/T  with  rcome  about  at  the  position  where  the  kft,  of  the  mode  jumps.  If  we 
could  have  assumed  a continuous  spectrum  in  k,  e<6/7ie  would  surely  be  a continuous  function 
of  r. 


We  now  discuss  what  effect  impurity  radiation  would  have  on  the  temperature  profiles. 
Experiments  on  TFR  have  shown  that  a large  fraction  of  the  input  ohmic  power  is  radiated 
away  by  impurities^®.  The  light  impurities,  typically  Oxygen,  radiate  away  the  largest  share  of 
the  power,  and  this  radiation  comes  from  the  edge  regions  of  the  plasma,  in  the  center  the  Ox- 
ygen is  fully  stripped and  there  can  be  no  line  radiation.  The  heavy  impurities,  typically 
molyiduim  radiate  away  a smaller  fraction  of  the  power  and  this  radiation  comes  from  the  cen- 
t.al  hot  region  of  the  plasma.  The  radiation  is  still  line  radiation;  molyiduim  is  not  fully 
stripped  until  the  plasma  reaches  a much  higher  temperature. 

We  will  first  discuss  the  impurities  near  the  center  where  there  is  anomalous  tranport. 
Our  basic  hypothesis  is  that  in  the  marginally  stable  regime,  the  plasma  picks  the  marginally 
stable  profile  nearly  independently  of  the  driving  mechanism.  Thus  if  there  is  an  additional 
energy  loss  term,  the  plasma  will  still  remain  at  a marginally  stable  temperature  profile,  but  the 
thermal  conduction  and  fluctuation  level  will  decrease.  However  if  the  number  of  impurities  is 
so  great  that  classical  conductivity  alone  is  sufficient  to  maintain  a stable  profile,  then  clearly 
this  profile  cannot  be  maintained.  If  the  number  of  impurity  ions  is  larger  than  this  critical 
value,  the  temperature  should  drop.  Thus  in  the  regions  where  there  is  anomalous  transport  a 
tokamak  electron  temperature  profile  should  be  relatively  insensitive  to  the  presence  of  impuri- 
ties as  long  as  the  number  of  impurities  is  sufficiently  small.  For  larger  concentrations  of  im- 
purities, the  electron  temperature  should  drop  sharply  as  the  impurity  concentration  increases. 

We  conclude  this  section  with  a discussion  of  the  effects  of  impurif-  radiation  for  plasmas 
whose  edge  region  is  classical.  If  the  edge  region  of  the  tokamak  is  classical,  and  the  electron 
and  ions  are  nearly  equilibrated  (T^,  — s T),  and  there  is  no  radiation  loss  or  energy  deposi- 
tion in  the  edge  region,  Eqs.  (I)  in  the  edge  region  can  be  reduced  to 
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r K,  ., 


II 

br 


constant. 


(21) 


Above  we  have  made  use  of  the  fact  that  A',.,  > > Since  the  ohmic  energy  channels  into 
the  hotter  regions,  it  is  reasonable  to  assume  no  energy  deposition  in  the  edge  region.  Equa- 
tion (21)  simply  states  that  the  total  energy  flux  through  each  cylinder  is  constant,  and  is  equal 
to  the  Ohmic  input  energy.  If  there  are  losses  due  to  radiation  (or  charge  exchange,  etc.),  then 
the  energy  flux  in  the  outer  region  must  be  smaller  because  some  of  the  ohmic  input  power 

has  been  radiated  away  before  it  reaches  the  edge.  This  then  implies  that  will  be  smaller 


at  each  temperature  in  the  classical  region,  and  thereby  lr(T)  will  decrease. 


br 


IV.  Results  for  Alcator 


Perhaps  the  most  striking  experimental  result  for  Alcator  is  the  increase  of  energy 
confinement  time  with  density'*^'*'*.  Our  code  does  predict  this  basic  dependence,  although  as 
with  TFR,  the  actual  confinement  time  is  too  large,  this  time  by  about  a factor  of  three.  In 
Fig.  9 is  shown  the  dependence  of  energy  confinement  time  on  central  density  for  three 
choices  of  field  and  current,  B — 50  KG,  / — 100  KA;  B — 75  KG,  / — 100  KA;  and  B — 75 
KG,  / • 150  KA.  In  all  cases,  Z^f  — 1 and  we  used  a hydrogen  plasma.  Notice  that  for  central 
densities  less  than  about  2 x lO''*  cm  the  confinement  time  increases  roughly  linearly 
with  density.  For  larger  tig  the  confinement  time  begins  to  decrease  again. 

Apparently  what  is  happening  is  that  as  long  as  a region  of  the  plasma  is  goverened  by 
anomalous  tranport,  the  central  temperature  and  'emperature  half  width  are  relatively  indepen- 
dent of  density.  Then  the  power  input  is  independent  of  density,  since  the  resisitivity  depends 
only  on  temperature.  The  energy  content  however  is  proportional  to  density,  so  that  the  ener- 
gy confinement  time  is  proportional  to  density.  When  the  density  is  sufficiently  high  however, 
the  plasma  is  in  the  Pfirsch  Schluter  regime,  and  there  is  no  trapped  particle  instability.  Then 
the  energy  confinement  time  begins  to  decrease  with  density.  For  all  three  choices  of  current 
and  field,  the  points  with  /?„  — 5 x fO*'*  were  completely  stable.  Not  only  are  there  no 
trapped  particle  instabilities  in  the  region  of  maximum  gradient,  q was  everywhere  larger  than 
one  so  that  are  no  MHD  modes  in  the  center  either.  Apparently  the  largest  confinement  time 
occurs  just  on  the  transition  from  anomalous  to  classical  thermal  conduction.  The  highest  elec- 
tron temperature  however  are  in  the  regime  of  anomalous  thermal  conduction  and  low 
confinement  time.  In  Fig.  10  a and  b are  plotted  the  radial  dependence  of  T,.  and  f,  for 
different  densities  for  the  case  B — 75  KG,  / • 100  KA.  Notice  that  at  low  density,  where  the 
thermal  conduction  is  anomalous,  the  temperature  profiles  are  quite  peaked.  At  higher  density, 
where  the  thermal  conduction  is  classical,  the  temperature  profiles  are  quite  broad. 

V.  Results  for  ATC 

Since  fluctuation  spectra  were  measured  in  ATC^*,  we  have  run  our  code  for  an  ATC 
plasma  in  the  uncompressed  state.  The  parameters  are  B„  — 16  KG,  7—65  KA,  ••  1.4 
X 10 Zgf  —2.  A plot  of  electron  temperature,  ion  temperature  and  17  as  a function  of  radius 
is  shown  in  Fig.  11.  In  Fig.  12  a and  b are  shown  the  radial  dependence  of  Ap,  and  efb/T^  — 
bn/itg.  Typical  values  of  eth/Tg  10  are  consistent  with  values  measured  by  microwave 
scattering  in  Ref.  45. 

VI.  Conventional  and  Unconventional 
Approaches  to  Electron  Heating 

Since  Ohmic  heating  has  not  heated  the  electrons  in  a tokamak  plasma  to  more  than 
about  2.5  Kev,  it  is  generally  agreed  upon  that  some  sort  of  supplementary  heating  is  necessary 
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lo  reach  reactor  temperatures.  Currently,  people  are  investigating  heating  by  microwave  and 
neutural  beams.  The  obvious  hope  is  that  by,  say,  doubling  the  input  power  one  can  double 
(or  more  than  double)  the  plasma  temperature.  For  our  purposes  here,  we  define  as  conven- 
tional heating  such  schemes  which  rely  on  increasing  the  input  power. 

However  there  are  two  sides  to  the  heating  problem,  the  energy  input  and  the  power 
losses,  ir  the  ion  losses  are  classical,  as  assumed  here,  then  there  is  no  way  they  can  be  re- 
duced. The  electron  losses  however  are  anomalous.  If  there  is  anyway  these  losses  can  be  re- 
duced, the  electrons  will  heat.  Furthermore  they  will  heat  without  additional  input  power,  so 
the  confinement  time  will  increase.  (Actually,  controlling  the  losses  probably  requires  input 
power  also).  We  will  define  as  unconventional  heating  such  schemes  which  rely  on  controlling 
the  electron  losses.  The  situation  is  rather  like  a phase  change*.  If  water  is  just  at  the  boiling 
point,  putting  in  excess  energy  does  not  raise  its  temperature.  Howevei  raising  the  boiling 
point  does  allow  the  water  to  be  heated  to  higher  temperatures.  The  problem  with  heating  by 
controlling  losses  is,  of  course,  that  there  is  now  no  general  agreement  on  what  causes  these 
losses,  much  less  on  how  these  losses  can  be  controlled.  Nevertheless,  at  least  our  model  of 
the  anomalous  losses  shows  that  there  is  a great  deal  to  be  gained  if  one  can  control  them. 

In  this  section  we  discuss  one  conventional  heating  scheme  and  three  unconventional 
heating  schemes.  All  runs  are  done  for  TFR  with  I — 140  KA  and  — 50  KG.  The  results 
are  summarized  in  Table  I. 

The  first  heating  scheme  we  examine  is  heating  of  the  plasma  by  neutural  beams.  We 
have  performed  two  runs.  The  first  of  these  runs  doubled  the  total  Ohmic  power  by  depositing 
38  KW  of  beam  power  into  the  ions  and  another  38  KW  into  the  electrons.  The  deposition 
profile  is  uniform  as  a function  of  radius.  As  is  shown  in  Table  1 the  central  electron  tempera- 
ture increases  from  about  2.5  Kev  to  about  2.8  Kev.  The  confinement  time  however  decreases 
from  108  ms  to  74  ms.  To  examine  the  effect  of  more  beam  power  and  a better  deposition 
profile,  we  have  done  another  run  tripling  the  original  Ohmic  power  to  76  KW  deposited  in 
ions  and  76  KW  in  electrons  and  choosing  an  energy  deposition  profile  proportional  to  the 
density  profile.  The  central  electron  temperature  now  rises  to  3 Kev,  but  the  confinement  time 
decreases  to  59  ms. 

One  can  understand  these  results  very  simply.  As  more  power  is  injected  into  the  plasma, 
the  temperature  gradient  tends  to  increase  and  the  plasma  tends  to  be  pushed  further  into  the 
plateau  and  banana  regime.  Thus  the  plasma  is  driven  into  a more  unstable  state.  However,  as 
the  plasma  is  pushed  into  a more  unstable  state,  it  is  fighting  the  thermostat  function,  which 
pushes  it  back  toward  stability.  The  harder  the  plasma  is  pushed  toward  instability,  the  harder 
it  fights  back.  Thus,  according  to  our  theory,  increasing  the  input  power  into  the  electrons  is 
not  a particularly  effective  way  to  heat  the  plasma. 

It  is  now  worth  while  to  compare  these  results  with  experimental  results  of  beam  injec- 
tion in  TFR  A quantitative  comparison  is  difficult  because  the  confinement  times  we  calcu- 
lated for  ohmically  heated  tokamaks  were  too  large  by  about  a factor  of  five.  Correspondingly, 
the  Ohmic  power  input  was  too  small  by  somewhat  less  than  this  factor.  These  are  some  quali- 
tative comparisons  that  can  be  made  however.  The  input  beam  power  into  TFR  was  also  com- 
parable to  or  larger  than  the  Ohmic  input  power,  and  the  power  deposition  was  about  half  into 
electrons  and  half  into  ions.  The  observed  ion  temperature  increased  by  about  I Kev.  Our 
code  showed  a considerably  smaller  increase  in  7",  (the  central  ion  temperature  only  increased 
by  about  200  ev).  The  reason  is  that  the  ions  are  in  the  Plateau  regime  where  the  thermal 
conductivity  is  an  increasing  function  of  temperature.  However,  once  the  ions  are  heated  to 
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the  bununa  regime,  where  the  thermal  conductivity  is  a decreasing  function  of  temperature, 
signiheant  heating  ion  should  result  with  increasing  input  power.  We  attribute  the  disagree- 
ment between  the  computed  and  experimental  ion  temperature  to  the  fact  in  the  experiment, 
much  more  power  is  going  into  the  ions  from  the  beams.  Thus  the  ions  can  be  heated  up  to 
the  banana  regime.  Additional  ion  heating  then  requires  much  less  beam  power. 

One  very  intesting  result  however  is  that  in  the  experiment,  the  electrons  either  do  not 
heat  at  all,  or  else  heat  by  no  more  than  10%.  Figures  8 and  13  of  Ref.  15  show  electron  tem- 
perature before  and  during  injection,  and  there  is  very  little  difference.  It  appears  that  the 
thermostat  function  controlling  the  experimental  temperature  profile  in  TFR  has  an  even  hard- 
er onset  than  the  thermostat  function  we  have  been  using,  Eq.  (7),  with  n —32. 

We  now  discuss  three  unconventional  heating  schemes.  First,  let  us  imagine  that  some 
way  could  be  found  to  stabilize  the  MHD  instability  at  the  q — 1 surface.  Since  this  is  a large 
scale  instability,  this  may  possibly  be  done  with  feedback  or  dynamic  stabilization.  It  is  in- 
teresting to  note  that  dynamic  stabilization  has  recently  been  shown  to  be  possible  for  the 
Raleigh-Taylor  instability  in  a laser  produced  plasma^.  While  we  do  not  know  if  dynamic  sta- 
bilization of  internal  kink  and  tearing  modes  is  possible,  our  results  show  that  it  is  worth  addi- 
tional investigation. 

In  our  numerical  simulation,  we  have  not  assumed  complete  stabilization  of  MHD  modes, 
rather  we  have  assumed  that  it  is  the  «/  — 1/2  surface  which  is  unstable  rather  than  the  q — 1 
surface.  The  result,  as  shown  in  table  1,  is  a very  impressive  amount  of  heating.  The  central 
electron  temperature  is  increased  to  4.4  Kev  and  the  confinement  time  increases  to  425  ms. 
Apparantly  what  is  happening  is  that  decreasing  t;(0)  by  a factor  of  two  increases  the  allowed 
shear  by  about  a factor  of  two.  (Recall  that  the  shear  is  related  to  the  change  in  l/q.)  The  in- 
creased shear  then  allows  for  larger  electron  temperature  gradients  before  reaching  the  thres- 
hold for  instability. 

Another  way  to  reduce  the  electron  losses  is  to  artificially  increase  the  collisionality  of  the 
electrons.  At  least  one  way  to  do  this  is  to  bathe  the  plasma  in  microwaves  at  the  electron  cy- 
clotron frequency"*’.  This  will  cause  enhanced  pitch  angle  scattering  of  the  electrons.  Thus 
the  temperature  at  which  the  electrons  enter  the  plateau  and  banana  regime  is  raised,  so  that  a 
larger  part  of  the  plasma  cross  section  will  be  stable.  Since  the  pitch  angle  scattering  is  caused 
only  by  the  microwave  power,  and  not  by  an  electron-ion  scattering,  the  plasma  resistivity  will 
be  unaffected,  so  there  will  be  no  additional  energy  input  caused  by  an  anomalous  resistivity. 

We  have  done  a numerical  simulation  in  which  the  collisionality  in  the  dispersion  rela- 
tion and  electron  thermal  conductivity  are  increased  by  a factor  of  2.5,  but  where  the  col- 
lisionality in  the  resisitivity  is  unchanged.  The  results  are  shown  in  Table  1.  The  central  elec- 
tron temperature  increases  to  2.9  Kev  and  the  energy  confinement  time  increases  to  145  ms. 

If  this  interpretation  is  correct,  it  implies  that  heating  experiments  at  the  electron  cyclo- 
tron frequency  will  do  as  well  or  better  by  injecting  power  with  a spectral  width  A oi/tu  * 2a/ R, 
instead  of  a single  frequency  chosen  to  match  the  cyclotron  frequency  at  a particular  point. 
Such  a spectral  width  would  allow  for  the  electron  collisionality  to  be  increased  across  the  en- 
tire plasma  cross  section. 

As  a final  example  of  the  control  of  electron  losses,  we  will  discuss  shaping  the  plasma 
cross  section.  It  has  recently  been  shown  that  an  elliptical  cross  section  both  slightly  increases 
the  shear  induced  damping"**  and  reduces  the  trapped  particle  induced  growth*^”**.  While  our 
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code  is  not  capable  of  handling  elliptical  cross  sections,  we  mock  up  the  effect  of  ellipticity  by 
increasing  the  (reduced.  See  Sect.  IlG)  shear  induced  damping  by  a factor  of  1.3  and  reducing 
the  trapped  particle  induced  growth  rate  by  a factor  of  two.  This  corresponds  roughly  to  an  el- 
lipiicity  of  three'***  As  is  shown  in  table  1.  the  central  electron  temperature  increases  to 
2.95  Kev  and  the  energy  confinement  time  increases  to  142  ms.  It  should  be  restated  now  that 
in  all  three  unconventional  heating  schemes,  the  power  used  to  reduce  the  electron  thermal 
conductivity  was  not  figured  in  the  definition  of  confinement  time.  The  reason  is  that  we  do 
not  know  what  this  power  is.  Surely  power  is  needed  to  stabilize  the  <y  — 1 surface  or  to  in- 
crease the  electron  collisionality  with  microwaves  at  the  electron  cyclotron  frequency.  Shaping 
the  cross  section,  however,  should  require  no  additional  power. 

To  summarize,  we  have  shown  that  conventional  heating  is  possible,  but  that  the 
confinement  time  (or  at  least  the  electron  confinement  time  in  experiments)  is  normally  re- 
duced. However  heating  by  controlling  the  anomalous  losses  offers  a potentially  very  attractive 
way  to  heat  the  electrons  in  tokamaks. 

VII.  Conclusions 

Our  conclusion  is  that  marginal  stability  theory  is  a very  useful  approach  to  the  problem 
of  electron  energy  confinement  in  tokamaks.  At  the  least  this  theory  provides  a worst  case 
analysis  of  the  anomalous  thermal  transport  resulting  from  a given  instability.  It  predicts  an 
energy  confinement  time  scaling  will  density,  as  seen  in  Alcator;  and  it  also  shows  that  the 
temperature  half  width  SriT)  and  the  radius  of  the  q — 1 surface  depends  only  on  the  single 
parameter  </(u).  as  seen  in  TFR.  Also,  it  provides  very  natural  explanations  for  the  failure  to 
heat  electrons  with  neutural  beams  in  TFR  and  also  for  the  heat  pulse  propagation  experi- 
ments in  ORMAK.  Finally,  it  suggests  that  controlling  the  electron  losses  may  be  a more 
effective  approach  to  electron  heating  than  injecting  additional  power. 
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APPENDIX 


The  functions  yj  and  defined  in  Sec.  II  contain  contributions  from  the  trapped  and  un- 
trapped electrons, 

/,  - (I  -i)  lm[T^  + K,|  (Al) 

/2  -{  Im{T2  + l^jl,  (A2) 

where 

K,  - r exp( 

' Jv  nV2  (1  -^11  R^)- 

y C exp(  -3/2) 

2 Jv  ^^2  (1  -^11  /j^) 

^0  Vf/ (Rt}<u)  -tu./w. 

j.  _ r exp  ( - 1^) 

* TT^'^  1 -I-  / /?,  f 

7 « f (f  ^ —3/2)  exp(  - 

^ 1 + //?,  ^ -3  f2 

Rx  ’‘i'J(elu)  “•  2^  e„  i>' (mif/Tm^)^^^lqb. 

U and  T denote  untrapped  and  trapped.  An  electron  is  considered  trapped  if  ^||  < (2  e) 
and  i > *'v^-  The  first  of  these  conditions  defines  the  boundary  of  the  loss  cone  in  velocity 
space  and  the  second  condition  results  from  demanding  that  the  effective  collision  frequency 
for  an  electron  is  less  than  the  bounce  frequency  in  the  magnetic  mirror. 

With  this  perscription  of  the  trapped  and  untrapped  regions  of  velocity  space  the  ima- 
ginary parts  of  the  untrapped  integrals,  f^x  and  can  be  obtained  by  straight  forward  integra- 
tion 

Im{  K,  1 — IT  exp  ( - Rg~^  ) [1  - exp  ( - v 1/^  ) !//?„ 

/„!  e*P  < ~ ) exp  ( - K ) ^ - /Rg 

where  we  have  assumed  Rg  > > t 

The  trapped  integrals  f)  and  can  be  further  simplified  by  performing  the  indicated  in- 
tegration in  spherical  coordinates,  /,  * 47r  (2e ) ; i^d(.  Thus,  the  functions  Tx  and 

‘■K* 

T2  contain  a transcendental  dependence  on  the  three  parameters  i>>g,  /?|,  and 

We  now  show  that  setting  the  lower  limit  in  the  T]  and  T2  integrals  to  zero  gives  a good 
approximation  to  these  integrals  for  all  If  ^ 1 so  that  the  lower  limit  of  integration  is 
significantly  different  from  zero  the  value  of  Rx  will  generally  be  quite  large  due  to  the  square 
root  of  the  mass  ratio  in  the  definition  of  X|.  Given  that  is  large  the  integrals  can  be 
evaluated  asymptotically  for  large  Rx . For  instance  the  asymptotic  value  of  Tx  is  found  to  be 


i(  ~ I’lf*  ) 1 -I-  »' + "T  '' V ■ 


NRL  MEMORANDUM  REPORT  3551 


The  explicit  dependence  of  (A3)  on  i .,,  is  actually  quite  weak  for  ~ 1.  If  >>  1,  in 
which  case  R\  >>  1,  and  T^  ■—  0,  even  if  the  lower  limit  of  integration  had  been  chosen  to 
be  zero  rather  than  <'],;■*,  both  expressions  are  so  small  that  either  one  gives  a good  approxima- 
tion to  the  necessary  shear.  Thus,  we  ignore  the  explicit  dependence  of  the  integrals  T^  and 
7’2  on  I’Y  by  setting  the  lower  limit  of  integration  to  zero. 


mials. 


The  remaining  dependences  of  T^  and  T2  on  /?,  and  were  approximated  by  polyno- 


17-,1 


, 1/2 


(2«) 


1/2 


£,  -R^ 


\ + a^  R^  + R}  + Rf 


L\h]  - 


, 1/2 


£2  + Ri  + \ .5R^ 

1 + fl2  ^1  + *2 


(A4) 


(A5) 


where  a^,  02,  6] , 62,  and  fi  are  functions  of  and  E\  and  £2  are  the  limiting  vaues  of 
/,„i7'| ) and  /'„|72'  as  £|  — 0 which  can  be  obtained  exactly  and  depend  on 

The  form  of  the  polynomials  and  the  dependences  of  aj,  b^,  02,  ^2’  and  0 on  ie„  were 
chosen  by  comparison  with  numerically  integrated  values  of  Ti  and  T2.  h can  be  seen  that 
the  approximations  (A4)  and  (AS)  agree  with  the  actual  values  of  T^  and  T2  for  both  large  and 
small  Afi . By  careful  choice  of  the  a’s  and  b's  we  were  able  to  get  agreement  within  10% 
between  the  polynomials  and  numerically  integrated  values  of  and  £2  for  the  relevant 
parameter  ranges  of  £|  and 


The  vaues  of  these  parameters  are 

- -1 
2 


'1 


1 

■y1 


and 


if  > 0.64  and 


exp  ( - 

(«€„) 

(A6a) 

exp(  - (ie„) 

(A6b) 

and  62  are 

1 -3.17 

- 7.264€„ 

(A7a) 

0 

0 

1 

(A7b) 

- 2.66  -1-  7.56 

(A7c) 

, - 11.18 

(A7d) 

■ -8.45  -1-  38.44 

(A7e) 

5.4  -2.65(^e„)/(l  -1-  570 ({«„ )■*). 

(A8a) 

8 ((«*„) 

- 0.2)  V ((««„)  +0.72). 

(A8b) 

+ (4.67 

+ 20.6(C«„))/(1  + 3.58  (4€„)^). 

(A8c) 

-0.98/(««„)  -1-  1.5, 

(A8d) 

• + 0.0025/ (C«„)2^ 

(A8e) 

for  < 0.64. 
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Table  1 


Scheme 

Added  Power 
Input 

r,(0) 

Only  Ohmic, 

P = 76  kw 

None 

2.5  kev 

108  ms 

Beam,  Uniform 
Deposition 

38  kw  into  ions, 
38  kw-electrons 

2.8  kev 

74  ms 

Beam, 

Deposition 

76  kw-ions 

76  kw-electrons 

3 kev 

59  ms 

Stabilization 
of  q = 1 

? 

4.4  kev 

425  ms 

Increasing 

Electron 

CoUisionality 

By  2.5 

2 

2.9  kev 

145  ms 

Cross  Section 
Shaping 

None 

2.95  kev 

142  ms 
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Fig.  1 — PLT  electron  temperature  profiles;  solid  line,  full  disper- 
sion relation  n = /iqCI  — -5  dashed  line,  full  dispersion  rela- 

tion n = «()( 1 — .85  r2/j2);  dash-dot  line  modified  dispersion 
relation  n = no(l  ~ -85  r^la^) 
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Fig.  3 — PLT  electron  temperature  profiles  obtained  with  (a)  classical 
and  (b)  neoclassical  resistivity 


Fig.  4 — PLT  electron  and  ion  temperature  profiles  showing  regions 
of  MHD  and  DTEM  instability 


22 


Fig.  5 — Temperature  half-width  and  width  of  the 
q{r)  = 1 region  in  TFR  as  a function  The 

solid  and  dash-dot  line  are  approximate  fits  to  the 
numerically  obtained  data  points.  The  dashed  line  is 
the  experimentally  observed  temperature  half-width. 


(a)  (b) 


’!»  r/« 


Fig.  6 - Electron  temperature  profiles  in  TFR.  (a)  numerically  obtained 
(b)  experimentally  obtained 
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Fig.  8 — (a)  The  most  unstable  mode  in  TFR  as  a 
function  of  r/n,  and  (b)  the  fluctuation  level  implied 
by  the  numerically  determined  anomalous  transport 
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